The Spin Stiffness and the Transverse Susceptibility of the 

Half-filled Hubbard Model 



Zhu-Pei Shi and Rajiv R. P. Singh 

Department of Physics, University of California, Davis, California 9561 6 

Abstract 

The T = spin stiffness p s and the transverse susceptibility x± °f the square 
lattice half-hlled Hubbard model are calculated as a function of the Hubbard 
parameter ratio U/t by series expansions around the Ising limit. We find that 
the calculated spin-stiffness, transverse susceptibility, and sublattice magne- 
tization for the Hubbard model smoothly approach the Heisenberg values for 
large U/t. The results are compared for different U/t with RPA and other nu- 
merical studies. The uniform susceptibility data indicate a crossover around 
U/t ~ 4 between weak coupling (spin density wave) behavior at small U and 

strong coupling ( Heisenberg ) behavior at large U. 
PACS numbers: 71.27.+a, 75.10.Jm, 75.40.Cx 
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Recent discovery of high-T c superconductivity in the Cuprate materials has generated 
tremendous interest in the subject of strongly correlated electrons. In these systems, the 
phenomena of unusual metallic behavior, antiferromagnetism and superconductivity occur 
in a narrow parameter range. It is widely believed that these phenomena have a common 
microscopic origin. The Hubbard model is one of the simplest models to describe correlated 
electron behavior in a solid. It consists of a single band of electrons, with nearest neighbor 
hopping parameter t and an on-site Coulomb repulsion between opposite spin electrons of 
magnitude U. The model is best understood at half-filling, that is, when there is one electron 
per unit cell, where the system becomes an antiferromagnetic insulator. At large values of 
U this system is well described by the Heisenberg model. At small U one expects the spin- 
density-wave (SDW) mean field description to become accurate. The possibility of d-wave 
superconductivity, away from half- filling, has been widely explored [|TJ . 

Direct calculations of the superconducting transition temperatures in the Hubbard model 
are beyond present numerical capabilities. Thus, phenomenological approaches, where the 
Hubbard model is used to determine the parameters in a scaling theory of superconductivity, 
are appropriate. If spin-fluctuations are important to the mechanism of superconductivity 
in the Cuprates, and if the magnetic excitations in the doped Cuprates are related to those 
in the stoichiometric insulating phases, the magnetic excitations in the half-filled system are 
clearly important for understanding superconductivity in these materials. 

The magnetic ground state of the two-dimensional (2D) square lattice Hubbard model has 
been investigated by quantum Monte Carlo simulations @-§J and Lanczos diagonalization || . 
These studies have mainly consisted of a finite size scaling analysis of the ground state 
properties of the model. They confirm the existence of long-range antiferromagnetic order 
at T = 0. To our knowledge, they have not been used to calculate the spin-stiffness constant 
or the spin-wave velocity for this model. Accurate numerical calculations of these quantities 
exist for the Heisenberg model. The Spin-density wave theory combined with the random 
phase approximation (RPA) has been used by Schrieffer et al. || to calculate the spin wave 
velocity for the Hubbard model. This calculation should become exact as U — > 0. Somewhat 
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surprizingly, it was found that the result of this calculation is also accurate in the Heisenberg 
limit {U/t > 1). 

In this Letter we first derive an expression for the spin stiffness constant of the Hubbard 
model by applying a slow twist in the ordering direction. We then introduce a one-parameter 
family of Hamiltonians which interpolate between the half-filled Hubbard and Ising models. 
This allows us to develop series expansions for the spin-stiffness. In addition, we develop 
series expansions for several thermodynamic parameters of this model, such as the uniform 
susceptibility and the sublattice magnetization. The spin wave velocity is calculated from 
the hydrodynamic relation @ u s 2 = p s /x±- 

At large U our results extrapolate smoothly to the Heisenberg values. They also show 
good agreement with the spin-density wave theory at small U. The variation of the magnetic 
susceptibility with the Hubbard parameter ratio U / 1 shows a relatively well defined crossover 
between axi~P behavior at large U and x± decreasing with increase in U at small U. 
This crossover between the strong coupling ( Heisenberg model behavior) and weak coupling 
SDW behavior occurs at U/t ~ 4. 

The Hubbard model is defined by the lattice Hamiltonian 



H = -t J2 (°ta Cja + C+ C ia ) + Uj2 ( n H 
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where cf a and Cj CT are the creation and annihilation operators for electrons with a 
z— component of spin cr at lattice site i, and n i)CT = cf a Ci a . U is the on site repulsive 
interaction, \i the chemical potential, and t the nearest-neighbor hopping amplitude. 

If we rotate the ordering direction by an angle 6 along a given direction such as y axis, 
then the spin stiffness constant p s can be defined through the increase of the ground state 
energy: E g {6) = E g {6 = 0) + \p s 6 2 + 0{6 4 ). This rotation can be carried out by the 
following transformation applied to the fermion operators: 



( 



cos (j) sin (f) 
—sin 6 cosd> 



(2) 



3 



After rotation by a relative angle 9 ( that is letting <f> change by 6/2 ) between neighboring 
sites separated along y axis (y denotes unit distance in y direction), H in Eq.(l) becomes 
H = H + H dia + H para + 0(9 3 ), where 



te 
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The "diamagnetic" term H dm is already of order 9 2 so for the calculation of the energy 
to order 9 2 , it can be replaced by its expectation value in the ground state of the 9 = 
Hamiltonian. We get pf a — — | (Kinetic Energy) — — | [E g {9 — 0) — ^ (n — L)], where n is 
a band filling, and L the local moment defined as L =< (n^ — n^) 2 >. The contribution of 
the "paramagnetic" term H para to the ground state energy in order 9 2 can be obtained from 
the expression, jf s ara = 2^f \e=o, where E is the energy of the Hamiltonian H + H para . 

In order to calculated these quantities numerically, we introduce an Ising anisotropy into 
the Hubbard Hamiltonian: 

H x = -\t Yl ( c taCi* + ci + ^(».t-!)K-!) +j(i-a) X! > 

<i,j>,a i ■ ~ <i,j> 

(3) 

where of = (rij-f — n^) is the z component of the spin at site i, and J a parameter which 
can be tuned to improve the convergence of the extrapolations. The particle-hole symmetry 
ensures half-filling. For A = the atomic limit of the Hubbard model is highly degenerate, 
however the Ising term selects from these the Neel states as the two degenerate ground 
states. Futhermore, this term also introduces a gap in the spectrum at A = 0. For A = 1 
the Ising anisotropy goes to zero and the conventional Hubbard Hamiltonian is recovered. 
Ground state properties of the model for A ^ can be obtained by an expansion in powers 
of A. If the gap does not close before A = 1 as expected for this model, we can obtain 
properties of the Hubbard model by extrapolating the expansions to A = 1. In the strong 
coupling limit, the half-filled Hubbard model is equivalent to the Heisenberg model with 
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Hamiltonian — Jh J2<ij> (o"t ■ Vj — 1) with Jh = t 2 /U. The optimum value of the 
parameter J is found to be near J H ||. 

We calculate series coefficients for the ground state energy and the local moment to 
11th order, and p para to 9th order in A. Using the Pade analysis, we obtain spin stiffness 
p s = 0.186(15), 0.15(1), 0.077(7) and 0.039(5) for U = 1,4, 10 and 20, respectively. They are 
plotted in Fig. 1 as filled circles with errorbars. The dashed line is a guide to the eye. The 
errorbars represent the spread in the Pade estimates. At U — 20, our result (p s U = 0.78) 
can be compared to that (0.73) of the Heisenberg model |]TT |. Our results for p s as well as 



M' are somewhat higher than the known results for the Heisenberg models. We believe this 
reflects the fact that the series have not converged as well as for the Heisenberg model and 
hence the reduction in these quantities due to the zero-point spin-wave fluctuations is not 
fully accounted for. Still, the convergence (~ 5%) is quite reasonable. 

We also compare our results with the Hartree-Fock approximation for the spin stiffness: 

Ps = -2iv £ k £(k g^ fc % where e ( k ) = ~ 2t ( cos k * + cos K) and W = t e2 ( k ) + A T /2 - A 
can be obtained by solving the gap equation ^ X) k jp + ^2)i/2 — 1- The spin stiffness given 

by mean field solution of the gap equation is plotted in Fig. 1 as a solid line. One can see 

that this approximation overestimates the stiffness at large U. At U = 20, the mean field 

result is p s U = 0.98 compared to 0.73 for the Heisenberg model. 

We note that the spin stiffness (filled squares in Fig. 1) from the variational Monte 
Carlo method with a Gutzwiller-type wave function are even larger than values of the mean 
field solution 0. We believe that the large discrepancy is due to the missing spin flip 
processes in the Gutzwiller variational wave function used in 0. Their calculations only 
get contributions to the spin-stiffness beyond the Hartree-Fock result from the "diamagnetic" 
term, whereas the "paramagnetic" part of the spin stiffness which contains spin-flip processes 
does not get corrected [[K] . 

The transverse susceptibility of the Hubbard model can be defined in the usual way as 
4x± = -2^1 \h=o, where E is the ground state energy of the Hamiltoniam H \ — h J2i °f • 
The expansion coeffients of Ax± to 9th order in A are also calculated. 



The transverse susceptibility obtained from Pade approximants are 4x± = 
0.94,0.58,0.75,0.92 and 1.40 for U = 1,4,6,10 and 20, respectively. They are plotted in 
Fig. 2 as filled circles. A dashed line is a guide to the eye. At U = 20, 4x±/U = 0.07 is close 
to the large U limit value x± J' — 0.065 for the Heisenberg model [Tl[ where J' = 4t 2 /U. At 
U — 4, our result (4%± = 0.58) also agrees with the Monte Carlo Simulations (~ 0.53) by 
White et al. ||. At U — 10, our result (Ax± = 0.92) agrees with the Monte Carlo simulations 
of Moreo [Q (0.98). One can see in Fig. 2. that the variation of the magnetic susceptibility 
with U changes character around U = U c ~ 4. It suggestes that there is a crossover at U c 
in the behavior of the 2D Hubbard model. For U > U c , the magnetic susceptibility (4%jJ is 
roughly proportional to U or inversely proportional to the spin superexchange J. This can 
be interpreted to imply that the magnetic state for U > U c for 2D Hubbard model is close to 
the Heisenberg model. But, for U < U c , the magnetic susceptibility (4x_l) has qualitatively 
different behavior. It decreases as U increases. This can be interpreted as the weak coupling 
behavior of the SDW antiferromagnetic state. 

We also compare our result with the mean field solutions of uniform magnetic sus- 
ceptibility: X ± = w [(1 - A 2 Uh -Ull/h)- 1 - 1], where h = ± £k ^ and 
I" 2 = ~Jn #3 (£) ■ K is plotted in Fig. 2 as a solid line. Again, we see that the mean field 
result is qualitatively correct, but overestimates the quantity especially at large U. 

The spin wave velocity obtained from the relationship (p s = v 2 x±) are v s = 
0.89, 1.02, 0.58 and 0.34 for U = 1, 4, 10 and 20, respectively. They are plotted in 
Fig. 3 as filled circles with a solid line as a guide to the eye. At U = 20, our v s = 0.34 agrees 
well with the limiting value (v s = 1.18\/2J' = 0.33) of the Heisenberg model [|TT . 



We compare this result with the RPA solution of Schrieffer, Wen and Zhang ||: v s = 
[4t 2 (l - A 2 Uh) g] 1 / 2 , where I 3 = ± £k ^ffi ■ It is plotted in Fig. 3 as a dotted line. 
We see that the RPA result for the spin wave velocity are lower than our series expansion 
results. 

To summarize, in this paper we have derived an expression for the spin stiffness of the 
Hubbard model, and calculated it and other key long-wavelength parameters by systematic 
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series expansion methods. Comparison with the 2D S — 1/2 Heisenberg antiferromagnet and 
RPA results, suggests that our extrapolations are good at small and large U. The uniform 
magnetic susceptibility data suggests a crossover between SDW behavior at small U and 
Heisenberg behavior at large U around U/t « 4. 
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Figure Captions 



Fig. 1 The spin stiffness as a function of U. Filled circles are our calculated data, with 
a dashed line as a guide to the eye. The errorbars indicate spread of the Pade approximants. 
The solid line is the mean field result. The filled squares are results of variational Monte 
Carlo simulations 0. 

Fig. 2 Transverse susceptibility as a function of U. Filled circles are our calculated 
data, with a dashed line as a guide to the eye. A crossover between SDW and Heisenberg 
antiferromagnet is suggested around U ~ 4. 

Fig. 3 The spin wave velocity as a function of U. Filled circles are our calculated 
data, with a solid line as a guide to the eye. The dotted line and dashed line represent 
results of RPA and Heisenberg antiferromagnet respectively. 
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